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We develop a kinetic theory of systems with long-range interactions taking collective effects and 

spatial inhomogeneity into account. Starting from the Klimontovich equation and using a quasilinear 

approximation, we derive a Lenard-Balescu-type kinetic equation written in angle-action variables. 

^D . We confirm the result obtained by Heyvaerts [Mon. Not. R. Astron. Soc. 407, 355 (2010)] who 

Cn ' started from the Liouville equation and used the BBGKY hierarchy truncated at the level of the 

fj \ two-body distribution function. When collective effects are neglected, we recover the Landau-type 

^ . kinetic equation obtained in our previous papers [P.H. Chavanis, Physica A 377, 469 (2007); J. Stat. 

^\ • Mech., P05019 (2010)]. We also consider the relaxation of a test particle in a bath of field particles. 

Its stochastic motion is described by a Fokker-Planck equation written in angle-action variables. We 

^*0 , determine the diffusion tensor and the friction force by explicitly calculating the first and second 

order moments of the increment of action of the test particle from its equations of motion, taking 

collective effects into account. This generalizes the expressions obtained in our previous works. We 

' ^ ' , discuss the scaling with N of the relaxation time for the system as a whole and for a test particle 

f\ • in a bath. 
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I. INTRODUCTION 



Systems with long-range interactions are numerous in nature. Some examples include self-gravitating systems, 

two-dimensional vortices, non-neutral plasmas and the Hamiltonian Mean Field (HMF) model (l|. The equilibrium 

C^ I statistical mechanics of systems with long-range interactions is now relatively well understood [2|. One important 

feature of these systems is that they may exhibit ensemble inequivalence and negative specific heats first observed 

in astrophysics (see [3|-|5[ for reviews). The out-of-equilibrium properties of these systems are also very interesting 

'^ ' [4]. In particular, long-range interacting systems may be stuck in non-Boltzmannian quasi stationary states (QSS) 
C , that persist for a very long time before finally reaching the Boltzmann distribution. These QSSs are stable steady 
states of the Vlasov equation which governs the "coUisionless" evolution of the system. In order to understand the 
different timcscales involved in the dynamics, it is important to develop a kinetic theory of systems with long-range 
interactions. Different kinetic theories have been proposed over the years. 

\N ■ Boltzmann Q developed a kinetic theory for dilute gases. In that case, the particles do not interact except during 

hard binary collisions. Boltzmann derived a kinetic equation from which he proved that entropy increases (iJ-theorem) 
^ and that the system ultimately relaxes towards the Maxwellian distribution of statistical equilibrium. 

^h , Landau [7| derived a kinetic equation for a neutral Coulombian plasma in which the charges interact via a 1/r 

potential. He started from the Boltzmann equation and considered a weak deflection approximation. Since the 
dynamical evolution of the system is dominated by weak binary collisions, it is possible to expand the Boltzmann 

f~^ ■ equation in terms of a small deflection parameter A ^ 1 and make a linear trajectory approximation. However, this 
treatment yields a logarithmic divergence at small and large impact parameters. The divergence at small scales is 
due to the neglect of strong collisions, resulting in large deflections. It can be cured by introducing a cut-off at the 
Landau length, corresponding to the impact parameter leading to a deflection at 90°. On the other hand, in a neutral 
plasma, the potential is screened on a distance corresponding to the Debye length [8[ . Phenomenologically, the Debye 

^^ • length provides an upper cut-off. Later on, Lenard [^ and Balescu [lOj developed a more precise kinetic theory that 
H I takes collective effects into account. This amounts to replacing the bare potential of interaction u(k) by a "dressed" 

- - ' potential of interaction w(k)/|e(k, k • v)| where e(k, k • v) is the dielectric function. Physically, this means that the 
particles are "dressed" by their polarization cloud. The original Landau equation, which ignores collective effects, is 
recovered from the Lenard-Balescu equation by setting |e(k, k- v)| = 1. However, when collective effects are taken into 
account, it is found that the logarithmic divergence that appears at large scales in the Landau equation is removed 
and that the Debye length is indeed the natural upper lengthscale to consider. At about the same time, Hubbard [ll| 
developed a test particle approach and derived a Fokker-Planck equation describing the relaxation of a test particle 
in a bath of field particles. He calculated the diffusion and friction coefficients by evaluating the first and second 
moments of the velocity dcfiection and took collective effects into account. The Fokker-Planck equation of Hubbard is 
recovered from the Lenard-Balescu equation in the "bath" approximation. These kinetic theories lead to a relaxation 
time towards the Maxwcll-Boltzmann distribution of the form tji ^ {A/lnA)t]j, where tjj is the dynamical time (the 
inverse of the plasma pulsation) and A the plasma parameter (the number of charges in the Debye sphere) . 



Chandrasckhar [12| - |14| developed a kinetic theory of stellar systems in order to determine the timescale of collisional 
relaxation and the rate of escape of stars from globular clusters. To simplify the kinetic theory, he considered an 
infinite homogeneous system. He started from the general Fokker-Planck equation and determined the diffusion 
coefficient and the friction force (first and second moments of the velocity increment) by considering the mean 
effect of a succession of two-body encounters. Since his approach can take large deflections into account, there is 
no divergence at small impact parameters and the gravitational analogue of the Landau length naturally appears 
in the treatment of Chandrasckhar. However, his approach leads to a logarithmic divergence at large scales that 
is more difficult to remove in stellar dynamics than in plasma physics because of the absence of Debye shielding 
for the gravitational force. In a series of papers, Chandrasckhar and von Neumann |15| developed a completely 
stochastic formalism of gravitational fluctuations and showed that the fluctuations of the gravitational force are given 
by the Holtzmark distribution (a particular Levy law) in which the nearest neighbor plays a prevalent role. From 
these results, they argued that the logarithmic divergence has to be cut-off at the interparticle distance. However, 
since the interparticle distance is smaller than the Debye length, the same arguments should also apply in plasma 
physics, which is not the case. Therefore, the conclusions of Chandrasekhar and von Neumann are usually taken 
with circumspection. In particular, Cohen et al. [l6[ argue that the logarithmic divergence should be cut-off at the 
Jeans length [ij], giving an estimate of the system's size, since the Jeans length is the gravitational analogue of 
the Debye length. Indeed, while in neutral plasmas the effective interaction distance is limited to the Debye length, 
in a self-gravitating system, the distance between interacting particles is only limited by the systems's size. These 
kinetic theories lead to a relaxation time of the form t/j ~ {N / In N)to where to is the dynamical time and N 
the number of stars in the system. Chandrasekhar [l^l also developed a Brownian theory of stellar dynamics and 
showed that, from a qualitative point of view, the results of kinetic theory can be understood very simply in that 
framework. In particular, he showed that a dynamical friction is necessary to reproduce the MaxwcU-Boltzmann 
distribution of statistical equilibrium and that the coefficients of friction and diffusion are related to each other by 
an Einstein relation (fluctuation-dissipation theorem). This relation is confirmed by his more precise kinetic theory 
based on two-body encounters. It is important to emphasize, however, that Chandrasekhar did not derive the kinetic 
equation for the evolution of the system as a whole. Indeed, he considered the Brownian motion of a test star in 
a fixed distribution of field stars (bath) and derived the correspondin g F okker-Planck equation [96| . This equation 
has been used by Chandrasekhar [1J|, Spitzer and Harm [20[, Michie [2l|, King [221, ^^'^ Lemou and Chavanis [2j] 
to study the evaporation of stars from globular clusters. King [24| noted that, if we were to describe the dynamical 
evolution of the cluster as a whole, the distribution of the field particles must evolve in time in a self-consistent 
manner so that the kinetic equation must be an integrodifferential equation. The kinetic equation obtained by King is 
equivalent to the Landau equation. There is, however, an important difference between stellar dynamics and plasma 
physics. Neutral plasmas are spatially homogeneous due to Debye shielding. By contrast, stellar systems are spatially 
inhomogeneous. The above-mentioned kinetic theories developed for an infinite homogeneous system can be applied 
to an inhomogeneous system only if we make a local approximation. In that case, the collision term is calculated as 
if the system were spatially homogeneous or as if the collisions could be treated as local. Then, the effect of spatial 
inhomogeneity is only retained in the advective (Vlasov) term which describes the evolution of the system due to mean- 
field effects [2^ [2g| . This leads to the Vlasov-Landau-Poisson system which is the standard kinetic equation of stellar 
dynamics. To our knowledge, this equation has been first written and studied by Henon [27|. Henon also exploited 
the timescale separation between the dynamical time tjj and the relaxation time tn 3> to to derive a simplified kinetic 
equation for /(e, t) , where e = w^/2 -I- <I>(r, t) is the individual energy of a star by unit of mass, called the orbit-averaged 
Fokker-Planck equation. In this approach, the distribution function /(r, v,i), averaged over a short timescale, is a 
steady state of the Vlasov equation of the form f{e,t) which slowly evolves in time, on a long timescale, due to the 
development of "collisions" (i.e. correlations caused by finite N effects or graininess). Cohn [2^ numerically solved 
the orbit-averaged Fokker-Planck equation to describe the collisional evolution of star clusters [93 • His treatment 
accounts both for the escape of high energy stars put forward by Spitzer [30[ , and for the phenomenon of core collapse 
resulting from the gravothermal catastrophe discovered by Antonov [31| and Lynden-Bell and Wood [sj on the basis 
of statistical mechanics. The local approximation, which is a crucial step in the kinetic theory, is supported by the 
stochastic approach of Chandrasekhar and von Neumann [l5| showing the preponderance of the nearest neighbor. 
However, this remains a simplifying assumption which is not easily controllable. In particular, as we have already 
indicated, the local approximation leads to a logarithmic divergence at large scales that is difficult to remove. This 
divergence would not have occurred if full account of spatial inhomo gen eity had been given since the start. The effect 
of spatial inhomogeneity was investigated by Severne and Haggerty [SJ] , Parisot and Severne [3J| , Kandrup [35| , and 
Chavanis [3a, |33]. In particular, Kandrup [33] derived a gene ralized Landau equation from the Liouville equation 
by using projection operator technics. Recently, Chavanis [36l . |37[ obtained this equation in a simpler manner from 
the BBGKY hierarchy or from a quasilinear theory, as an expansion in 1/A^ in a proper thermodynamic limit. This 
generalized Landau equation is interesting because it takes into account effects of spatial inhomogeneity and memory 
which were neglected in previous approaches. Since the finite extension of the system is properly accounted for, there 



is no divergence at large scales [gg. Furthermore, this approach clearly shows which approximations are needed in 
order to recover the traditional Landau equation. Unfortunately, the generalized Landau equation remains extremely 
complicated for practical applications. In addition, this equation is still approximate as it neglects "collective effects" 
and considers binary collisions between naked particles. As in any weakly coupled system, the particles engaged in 
collisions are dressed by the polarization clouds caused by their own influence on other particles. Collisions between 
dressed particles have quantitatively different outcomes than collisions between naked ones. In the case of plasmas, 
collective effects are responsible for Debye shielding. For self-gravitating systems, they lead to "anti-shielding" and are 
more difficult to analyze [931 ■ Some authors like Thorne [33 , Miller [40] , Gilbert [4l|, \^ , Lerche [4^] , and Weinberg 
|44| attempted to take both spatial inhomogeneity and collective effects into account. However, they obtained very 
complicated kinetic equations that have not found application until now. They managed, however, to show that 
collective effects are equivalent to increasing the effective mass of the particles, hence diminishing the relaxation time. 
Since, on the other hand, the effect of spatial inhomogeneity is to increase the relaxation time [SJ], the two effects act 
in opposite direction and may balance each other. 

The same difficulties (and addition al on es) occur for other systems with long-range interactions [2[. For example, 
the kinetic theory of the HMF model [lOOl | is v ery complicated. In the homogeneous phase (at supercritical energies) 
the Lenard-Balescu collision term vanishes [43, 14^, like in the case of homogeneous ID plasmas J4^|50l, ^^ the kinetic 
theory must take higher order correlations (three-body, four-body...) into account. In the inhomogeneous phase (at 
subcritical ene rgies), we have to face the problems of spatial delocalization discussed above. The kinetic theory of 2D 
point vortices J5lM53l | presents very similar features. Finally, ID stellar systems are always spatially inhomogeneous 
and it is not possible to make a local approximation, contrary to the 3D case. 

An important problem in the kinetic theory of systems with long-range interactions is therefore to derive a kinetic 
equation that takes both spatial inhomogeneity and collective effects into account. A standard tool to deal with 
spatially inhomogeneous systems is to use angle-action variables. This has been used by many authors in astrophysics 
in order to solve dynamical stability problems l54l458l | or to compute diffusion and friction coefficients in a Fokkcr- 
Planck approach [59l - l66| . In an earlier paper [63|, starting from the Klimontovich equation and using a quasilinear 
approximation, we attempted to derive a Lenard-Balescu-type equation with angle-action variables, taking collective 
effects into account. However, in the course of the derivation, we made a factorization assumption that was not 
justified. We had the intuition that this problem might be solved by expanding the density-potential functions on a 
"good" basis but we did not succeed in developing an adequate formalism. However, we argued that the final kinetic 
equation that we obtained (sec Eq. (37) of (63) is correct when collective effects are neglected. In a subsequent paper 
|68| , we derived this kinetic equation in a more satisfactory manner without ad hoc assumption. Since collective effects 
are neglected, this kinetic equation can be viewed as a generalization of the Landau equat ion, written in angle-action 
variables, applying to spatially inhomogeneous systems with long-range interactions [lOl| . 

In the meantime, Heyvaerts [69| managed to derive a Lenard-Balescu-type kinetic equation with angle-action 
variables starting from the Liouville equation and using the BBGKY hierarchy truncated at the level of the two- 
body distribution function (i.e. neglecting three-body correlations). To calculate the collective response, he used 
a biorthogonal basis of pairs of density-potential functions. This is a relatively classical tool in astrophysics |7C 



introduced by Kalnajs [7l|, that we did not know at the time of our papers [67|,l68[. Although the derivation given 
by Heyvaerts [63 is very satisfactory, we will show in the first part of the present paper (Sect. HBt hat the Lenard- 
Balescu-type kinetic equation can also be derived from the Klimontovich approach initiated in [67| . The use of a 
biorthogonal basis of pairs of density-potential functions solves the shortcomings of our previous derivation. This new 
derivation may be technically simpler than the one given by Heyvaerts [69[ . Furthermore, this derivation is completely 
parallel to the one exposed in [2, UA, UM to derive the ordinary Lenard-Balescu equation and in [5l|, [Sj] to derive 
a Lenard-Balescu-type equation appropriate to two-dimensional point vortices. Our approach provides therefore a 
unified description of kinetic theories for systems with long-range interactions when collective effects are accounted 
for. This generalizes the kinetic theory presented in [68| when collective effects are neglected. Our quasilinear theory 
has some similarities with the theories developed by Weinberg |6J| , Ma and Bertschinger [TJ] : and Pichon and Aubert 
j6a | in astrophysics. However, in our case, the source of noise is due to discreteness (finite N) effects internal to the 
system while in [64l . l65l . |74| | it is caused by external sources (perturbations on a galaxy, cosmological environment on 
dark matter halos...). In the second part of the paper (Sect. IIIip . we consider the relaxation of a test particle in a 
bath of field particles. Its stochastic motion is described by a Fokker-Planck equation with angle-action variables. 
We determine the diffusion tensor and the friction force by explicitly calculating the first and second order moments 
of the increment of action of the test particle from its equations of motion, taking collective effects into account. 
This generalizes the expressions obtained in our previous works [GTI. l68j. neglecting collective effects. We show the 
connection between the Fokker-Planck equation describing the relaxation of a test particle in a bath and the kinetic 
equation describing the evolution of the system as a whole. We also discuss the scaling with N of the relaxation time 
for the system as a whole and for a test particle in a bath. 



II. COLLISIONAL EVOLUTION OF THE SYSTEM AS A WHOLE 
A. The Klimontovich equation 

We consider an isolated system of material particles in interaction moving in a d-dimensional space. Their dynamics 
is described by the Hamilton equations 

dri dH dvi dV. ,., 



with the Hamiltonian 

H 



^ 1 
Y^-mv^ + m^Y,ui\r,-rj\), (2) 



2 

1=1 i<j 



where M(|r — r'|) is the potential of interaction and m the individual mass of the particles. For an isolated system, 
the total energy E = T-L is conserved. We assume that the potential of interaction decays at large distances like r~^ 
with 7 < d. In that case, the potential is said to be long-ranged Q]. Long-range potentials include the gravitational 
and the Coulombian potentials, corresponding to 7 = d — 2. 

We introduce the discrete distribution function /^(r, v,i) = m^^6{Y — Yi{t))5{-v — Vi(f)) and the corresponding 
potential $d(r,0 = /""(k ~ r'|)/rf(r', v',i)(ir'(iv'. Using the equations of motion ([IJ and ([5]), we find that /^(r, v,i) 
satisfies the Klimontovich [T^ equation 

^ + [Hd,/rf]=0, (3) 

where i/^ = v'^ /2 + $d(r, t) is the exact energy of an individual particle and we have introduced the Poisson bracket 

^^f dr dr dv 

The Klimontovich equation ^ contains exactly the same information as the Hamiltonian equations (HJ-Q. 

We now introduce a smooth distribution function /(r, v, t) ~ (/d(r, v, t)) corresponding to an average of /rf(r, v, t) 
over a large number of initial conditions. Wc then write /^ = f + Sf where Sf denotes the fluctuations around 
the smooth distribution. Similarly, we introduce $(r, f) = ($c;(r,t)) and write <f>d = $ + (5$. Substituting these 
decompositions in Eq. ([3]), we obtain 

^ + ^ + [i^, /] + [H, Sf] + [S^, /] + [S<i>, Sf] - 0, (5) 

where we have written Hd = H + SH with i/ = ti^/2 + $ and SH = (5$. Taking the ensemble average of Eq. (O, we 
get 

^^+[H,f]^-m,Sf]). (6) 

Subtracting this relation from Eq. ([5]) and neglecting the nonlinear terms, we obtain an equation for the perturbation 

^ + [H,Sf] + [S<^J]=Q. (7) 

This linearization corresponds to the quasilinear approximation J72| . It is valid at the order \/N in a proper thermo- 
dynamic limit A^ — > -|-(X) in which the system size scales like V r^ \, the individual mass like m r^ 1^ the total energy 
like E ^ N, the inverse temperature like /3 ~^ 1, the dynamical time like f/3 ~ 1, and the potential of interaction like 
u ^ \/N [7i|. This is equivalent to putting 1/N in front of the potential energy in the Hamiltonian (Kac prescrip- 
tion). In this thermodynamic limit, the distribution function and the potential are of order f /N ^ $ ^ 1 and their 
fluctuations are of order Sf /N ~ i5$ ~ 1/VN [331 ■ Therefore, the term in the right hand side of Eq. (|6]) is of order 
1/A^ (as compared to the left hand side). For N — > +00, it vanishes and we obtain the Vlasov equation 

^ + [HJ]^0. (8) 



In order to deal with spatially inhomogeneous systems, it is convenient to use angle-action variables [70|, |76[. By 
construction, the Hamiltonian H in angle and action variables depends only on the actions J = (Ji, ..., J^) that are 
constants of the motion (the conjugate coordinates w ~ {wi, ..., Wd) are called angles). Therefore, any distribution of 
the form / = /(J) is a steady state of the Vlasov equation. According to the Jeans theorem, this is not the general 
form of Vlasov steady states. However, if the potential is regular, for all practical purposes any time-independent 
solution of the Vlasov equation may be represented by a distribution of the form / = /(J) (strong Jeans theorem). We 
assume that the system has reached a quasi stationary state (QSS) described by a distribution of the form / ~ /(J) as 
a result of a violent collisionless relaxation [73| involving only meanfield effects. Due to finite N effects (granularities, 
correlations...), the distribution function / will slowly evolve in time. This is similar to the effect of collisions in the 
theory of gases. Finite TV effects are accounted for in the right hand side of Eq. © which is similar to the collision 
term in the Boltzmann equation. This term is expected to drive the system towards statistical equilibrium. Therefore, 
it is expected to select the Boltzmann distribution among all possible steady solutions of the Vlasov equation. Since 
this term is of order 1/A^ or smaller, the effect of "collisions" (granularities, finite N effects, correlations...) is a very 
slow process that takes place on a timescale Ntn or larger. Therefore, there is a timescale separation between the 
dynamical time tu ^ 1 which is the timescale at which the system reaches a steady state of the Vlasov equation 
through violent collisionless relaxation and the "collisional" time tn at which the system is expected to relax towards 
statistical equilibrium due to finite N effects. Because of this timescale separation, the distribution function will evolve 
through a sequence of stationary solutions of the Vlasov equation depending only on the actions J, slowly changing 
with time due to "collisions" (finite N effects). Hence, the distribution function averaged on a short (dynamical) 
timescale can be approximated by 

/(r,v,i)^/(J,t). (9) 

In that case, the distribution function is a function / = /(J,t) of the actions only that slowly evolves in time under 
the effect of "collisions" . This is similar to an adiabatic approximation. The system is approximately in mechanical 
equilibrium at each stage of the dynamics and the "collisions" slowly drive it to thermodynamical equilibrium. The 
purpose of this paper is to obtain a closed expression for the collision term in action space when the rapid dynamics 
has been averaged over the angles. 

B. Angle-action variables 

The previous decomposition has separated the discrete distribution function /^ into a smooth component / that 
slowly evolves with time and a fluctuating component Sf that changes more rapidly. When we focus on the evolution 
of the fluctuations, we can consider that the smooth distribution is "frozen" . The smooth distribution / determines a 
smooth potential $ that, we assume, is integrable. Therefore, to analyze the evolution of the fluctuations, it is possible 
to use angle-action variables constructed with this smooth potential [70, IZyl • A particle with coordinates (r, v) in 
phase space is described equivalently by the angle-action variables (w, J). The smooth Hamiltonian equations for the 
conjugate variables (r, v) are 



rfr _dH dv _ dH 

dt dv ' dt dr 



(10) 



If we use the variables (r,v), wc find that the dynamics is relatively complicated because the potential explicitly 
appears in the second equation. Therefore, this equation dv/dt — — V$ cannot be easily integrated except if $ = 0, 
i.e. for a spatially homogeneous system. In that case, the velocity v is constant and the unperturbed equations of 
motion reduce to r = vi -I- Tq, i.e. a rectilinear motion at constant velocity. Now, the angle-action variables are 
constructed so that the Hamiltonian does not depend on the angles w. Therefore, the smooth Hamiltonian equations 
for the conjugate variables (w, J) are 

dw dH ^,,, dJ dH 

where n{3) is the angular frequency of the orbit with action J. From these equations, we find that J is constant and 
that w = n{3)t -\- wq. Therefore, the equations of motion are very simple in these variables. They extend naturally 
the trajectories at constant velocity for spatially homogeneous systems. This is why this choice of variables is relevant 
to develop the kinetic theory. Of course, even if the description of the motion becomes simple in these variables, 
the complexity of the problem has not completely disappeared. It is now embodied in the relation between position 
and momentum variables and angle and action variables which can be quite complicated. In this paper, we shall not 
consider this problem and we shall remain in action space. 



Since the Poisson bracket is invariant on a change of canonical variables, we have 



83 dw dw 83 
Now, in terms of the angle-action variables, using the relations 



dH 8f 8H 8f 
[■«, Jj - ITT • irr " 7^ • 777- U^) 



8H 85 f 8H 85 f 85 f 

[..,/]^^-|^-^-|^-f5-|, (14) 

83 aw aw 83 aw 83 



8f 8 /.rd5^\ ,,,. 



and averaging over the angles, we find that Eqs. (H])-© take the form 

8t 83 \ dw 

These equations govern the evolution of the smooth distribution function and of the fluctuations in angle and action 
space at the order 1/A^. 

C. Lenard-Balescu and Landau-type kinetic equations 

Equations (|15p and ((T6)) form the basis of the quasilinear approximation. We shall assume that the fluctuations 
evolve rapidly compared to the transport time scale, so that time variation of f{3,t) and J7(J,i) can be neglected in 
the calculation of the collision term (Bogoliubov ansatz). In other words, the distribution function can be regarded 
as constant in Eq. (|16p because it evolves on a (relaxation) timescale which is much longer than the (dynamical) time 
required for the correlation function {5f5^) to reach an equilibrium for a given value of /. This adiabatic hypothesis 
is valid for N ^ 1. For brevity, we shall omit the variable t in / and O. It will be restored at the end, when we 
take into account the time dependence of the smooth distribution function through the kinetic equation (jlSp . We 
also assume that the distribution function / remains Vlasov stable, so it evolves under the sole effect of "collisions" 
and not because of dynamical instabilities. In that case, Eqs. (fT5|) - (fT6)) can be solved with the aid of Fourier-Laplace 
transforms and the collision term can be explicitly calculated. We present below the calculations that lead to the 
Lenard-Balescu-type equation (|56p for spatially inhomogeneous systems. 

The Fourier-Laplace transform of the fluctuations of the distribution function 5f is defined by J102J : 

J^l rfte-^''"'-"*''5/(w,J,i). (17) 

This expression for the Laplace transform is valid for Im(w) sufficiently large. For the remaining part of the complex 
Lj plane, it is defined by an analytic continuation. The inverse transform is 



<5/(w, 3,t)^Y.lT- e'^''-"-"*^'^/(k, J, ^), 
^ Jc 27r 



(18) 



where the Laplace contour C in the complex ui plane must pass above all poles of the integrand. Similar expressions 
hold for the fluctuations of the potential (5$(w, J,t). If we take the Fourier-Laplace transform of Eq. ([TC|. we find 

that 

8f ~ 
- Sf(k, J, 0) - ioj (5/(k, 3,Lu)+ik-n 5/(k, J, w) - ik ~ 5$(k, J, w) = 0, (19) 

where the first term is the spatial Fourier transform of the initial value 

SfiK J' 0) - / ^ e-*-5/(w, J, 0). (20) 



The foregoing equation can be rewritten 



,;>,j.„).jLiL,4„.a,.). Wk.J,o, 



kJ7- 



LJ 



z(k ■ ft — Lj) 



(21) 



where the first term on the right hand side corresponds to "collective effects" and the second term is related to the 
initial condition. In order to simplify the notations, we have noted O for 0(J). Similarly, we shall note / for /(J), 
ft' for ri(J') and /' for /(J')- The fluctuations of the potential $(r, t) are related to the fluctuations of the density 
p(r, i) = //(r, v,i)dv by 5$(r, t) = J u{\r — r'\)6p{r' ,t) dr' . We introduce a bi-orthogonal basis [701 i^i which the 
density Pa{r) and the potential $a(r) are related to each other by 



$a(r) ^ fu{\r - r'|K(r')dr', /"p,(r)$;,(r) dr = -J„,,.. 

The fluctuations of density 6p{r,t) and potential (5$(r,i) can be expanded on this basis as 

dp{r, t) = Y, Mt)Pa{r), <5$(r, t) = ^ A^{t)<i>^{r). 

a a 

Introducing angle-action variables, we have 

5$(w,J,t) = ^yl„(t)$„(w,J). 

a 

Taking the Fourier-Laplace transform of this equation, we obtain 

,5$(k,J,w) = ^i„(u;)$„(k,J), 



where 



p + oo 

i„(w)= / A^{t)e"^Ut, 
Jo 



is the Laplace transform of Aa{t) and 

$„(k,J) = y$„(w,J)e-'''-- 
is the Fourier transform of $q(w, J). We can rewrite Eq. ([2T|) as 

Sf{k,3,uj 



k-^ 

— ^f ■<- — ^ 



dw 

(2^' 

5/(k,J,0) 



k. ■ ft — uj 

a 

Taking the inverse Fourier transform of this expression, we obtain 

k.^ 



z(k ■ ft — ijj) 



Sf{w,3,u;) = Y, e^''--<5/(k, J,c.) = ^ e*- 



9J 



]i ■ ft CJ 



^i„H$„(k,j) 



k k 

Multiplying the left hand side of Eq. ([^^ by $*,(w, J) and integrating over w and J, we get 



'^/(k,J,0) 

i(k ■ ft — uj) 



(5/(w,J,a;)$;,(w,J)dwdJ= / 5/(r, v, w)$;,(r) drdv = / Jp(r, w)$;,(r) dr 

^H / ^a(w)pa(r)$;,(r)dr = -^i„(w)(5a,a' =-ia'(w), 



(22) 



(23) 



(24) 



(25) 



(26) 



(27) 



(28) 



(29) 



(30) 



where we have used the fact that (w, J) are canonical variables so that drdv = dwdJ. Multiplying the right hand 
side of Eq. ([^ by $*, (w, J) and integrating over w and J, we get 



5:/dwdJe--$;,(w,J) 
k -^ 

= (2,)''j;/<iji;,(k,j) 



''^a j;i„(.,4„(k,j, , «>'.■'■«) 



k. ■ ft — ui 

k.^ 



dJ 



k-ft 



^i„(c^)d,(k,j) 



i(k • 17 — oj) 

Sf{k,3,0) " 
j(k -17 — 0;) 



(31) 



From Eqs. ^Q and (|3Tj), we obtain 



i„,H = (27r)'^^ 



1. ^ 



a L k 

Introducing the dielectric tensor 



j-j-ci.(k,j)4:,(k,j) 



i„M + (2^)^^|dJ^ffll^<|.:,(k,J). (32) 



eaa'(w) = -^oa' + i^nf J2 / ^J j^ . ^^:' ^ ^^ (k, J)4„, (k, J) , 



the foregoing relation can be rewritten |l03l | : 



Inverting the matrix e(a;), we get 



'5/(k,J,0) 



[k-ft-uj] 



<5:(k,j). 



i„M ^ -i2nrY.(^-XMY^ J dj' ,yy:°> :'(^^^ J' 



Finally, substituting Eq. ([35]) in Eq. psp and introducing the notation 

1 



i?k,k,(J,J',Lj) 



^|.„(k,J)(e-iWH$*,(k',J'), 



we obtain 



.<.(k,J,.)^-(2.)^E/^J'l^^^ 



_^/(kViM3^ 

w) i(k' • fJ' - w) 



(33) 



(34) 



(35) 



(36) 



(37) 



We can use the previous expressions to compute the collision term appearing on the right hand side of Eq. ([15 
One has 



Sf 



i9w 



did 



= y ^y i!^,k' 



ir 27r ^-—^ Ir 27r 



dio' 



■i(k-w— cjt)i(k'-\v — cj'i) 



{Sf{-k,J,u;)S^k',3,iu')). 



Using Eq. (|2T|), we find that 

{Sf{k,J,u;)S^k',J,u;')) 



^ • - -(^^(k, J, .)^$(k', J, .')> + ^'^"^^' '' "^'*^^' '' "'^^ 



k • J7 — oj 



i(k • Jl — a;) 



(38) 



(39) 



The first term corresponds to the self-correlation of the potential, while the second term corresponds to the correlations 
between the fluctuations of the potential and of the distribution function at time t ~ 0. Let us consider these two 
terms separately. 

From Eq. (|37p . we obtain 



(5$(k,J,a;)jl>(k',J,a;')) = -(27r)2'^ ^ f dJidJ 

ki,k2 



1 



1 



((5/(ki,Ji,0)(5/(k2,J2,0)) 



i?k,ki(J, Jl,w) Dk'.k.(J, J2,C^') (kl • f2i - W)(k2 • ^2 - LO') 



Using the expression of the auto-correlation of the fluctuations at t = given by (see Appendix [C|) : 



(<5/(k,J,0)(5/(k',J',0)) = 4,_k"5(J - J')W(J), 

(ztt)" 



we find that 



{6^{k,3,u;)d^k',J,u;')) = (2^)^^ ^ f dJ^ ^ 



m/(Ji 



w) Z?k',-ki(J, Ji,w') (kl • f^i - w)(ki • fii + w') 



(40) 



(41) 



(42) 



Considering only the contributions that do not decay in time, it can be shown [T^] that [(ki • Jli — w)(ki • fti + uj')] ^ 
can be substituted by {2tt)'^S{uj + uj')S(ki ■ fti — u). Therefore 

{5nk,3,co)5^k',3,Lu')) ^ {2nf+^mS{c, + co')y] fdJ,- ^1- j\- -S{k^-n^-u)f{3,). (43) 

Similarly, one finds that the second term on the right hand side of Eq. f|39p is given by 

(^/(k,J,0)Jl»(k^J,c.O) 1 o ^f(^^ (aa^ 
-— — = - 27r m5{uj + uj)- — -dk-O-tj/ J). (44 



From Eq. (|43)) . we get the contribution to Eq. (jSSj) of the first term of Eq. (|39p . Since the colhsion term depends 
only on the action, it can be averaged over w without loss of information. This br ings a Kronecker factor i5k,-k' which 
amounts to taking k' = — k. Then, using Z3_k.-ki(J, Ji, —^) = ^k.ki(J, Ji, w)* [iQJ], we obtain 

^/^\ ^^^ii-Kf+'mY I — I di'k-^-^, J__^5(k' . f2' - c.)/(J'). (45) 

aw/, ^ ' ^Jc^ttJ k-rj-^|i5k.k'(J,J',^)P ^ ^^^ ' ^ ' 

Finally, using the Landau prescription w — > oj + tO+ and the Plemelj formula 

^ v(-]TiTrS{x), (46) 



a; ± iO+ \x 

where V denotes the principal value, we can replace l/(k • fi — cj — iO^) by +iTrS{'k ■ ft — lu). Then, integrating over 
Lo, we obtain 

From Eq. (|44p . we get the contribution to Eq. (|38| of the second term of Eq. (p9|) . It is given by 



Sf- 



)^rV''H D..i3,j,k.n) )-^^'^- (^«) 



d^^ii k 

Using Eq. p6p and introducing the convention of summation over repeated indices, we have 

1 



2?k,k(J,J,krj) 

Therefore 



$„(k,J)(6-^)„„,(k.n)$;,(k,J). (49) 



Im ( ^^^(j J j^^^ ) = ^*"(k, J) [{e-'U'{k ■ n) - {{e-'U,^{k ■ n))*] K'{K3). (50) 

To compute the term in brackets, we use the identity 

(e-i)_(e-i)t = ,-i(et_,)(,t)-i. (51) 

In component form, it can be written 

{e-%c.' - iie-%'c.r - ie-')ax{el,^ - eAA')((e"')a'A')*- (52) 

From the expression l\'3'3\i of the dielectric tensor, we find that 

eA'A(w)*-eAA'(w) = -^(27^)'^+l^ /dJ ("k • 1^") ,5(k • fJ - w)$^(k, J)$A'(k, J), (53) 



where we have used the Landau prescription uj ^- uj + iO+ and the Plemelj formula (j46]). Inserting Eqs. ([52]) and ([53 
in Eq. ([50]) and using Eq. (|36p . we obtain 



Ini ( -^ ) = -7r(27r)'* V / rfJ' ^ —S(k -U-k' -n') (k 

Vi?k,k(J,J,k-n); ^ ' ^J |i5k,k'(J,J',k-f2)P ^ '\ 



'il- («) 
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Substituting Eq. (HH) in Eq. dH]), wc find that 

(^^^)„ = --(2-)'- g / ''' "^ ii....(j,JMc.o)P ^('- • ^ - '^' ■ ^') (^' • i) ^(•^)- (^^) 

Finally, regrouping Eqs. (|15p . (j47p and (j55p and restoring the time variable, we end up with the kinetic equation 

k,k' I : \ ^ y I \ / 

This Lenard-Balescu-type kinetic equation, taking spatial inhomogeneity and collective effects into account, was 
derived by Heyvaerts [63] from the BBGKY hierarchy. It describes the "collisional" evolution of the distribution 
function caused by the weak noise created by the discreteness of the particles accompanied by their associated 
polarization cloud. It is written in action space, which is possible when the Hamiltonian associated with the average 
potential is integrable. For spatially homogeneous systems, it reduces to the ordinary Lcnard-Balescu equation 

I - <'-r-^ • / '^^^^' "^ Ki^'t' • ^^ - ^'^^' il-^) ^(^' ^)^(^' *) (^^) 

where 

e(k, w) = 1 + (27r)'^M(k) / ^ dv, (58) 

J uj — k • V 

is the dielectric function. When collective effects are neglected, we obtained in Refs. [63, |6^ a Landau-type kinetic 
equation with angle-action variables of the form 

^ = Ti{2i,Ym^ ■ Y,ldJ'k\Au,u'{J,J')\^S{k . O - k' . O') (^k ~ _ k' . A^ /(J,t)/(J',t), (59) 



k,k 



where Ak.k'(J, J') is the Fourier transform of the potential of interaction written in angle- action variables (see Ap- 
pendix [X|. This equation was derived in [6q by a method that does not use Fourier-Laplace transforms nor biorthog- 
onal basis. It arises from the generalized Landau equation 



df df N-1 

dt dr N 



^^■l^ = i^l/^l'^'^'^^^^'^'^^^^('^'-^^ 



■^^a->o)A+^^(o_i)^ 



/(r,v,t-r)^(ri,vi,t-r), (60) 



derived in 35i-l37ll . which is valid for systems that are not necessarily spatially homogeneous and not necessarily 
Markovian [l05l | . The Landau-typ e eq uation (|59p is obtained from Eq. (|60p by making a Markovian approximation and 
introducing angle-action variables |68| . For spatially homogeneous systems, it reduces to the ordinary Landau equation 
corresponding to Eq. ([57|) with |e(k, k • v)| = 1. The connection between Eqs. (|56p and (|5^ is clear. If we neglect 
collective effects, l/_Dk.k'(J, J',w)bare is just the Fourier transform Ak,k'(J' J') of the "bare" potential of interaction 
(see Appendix \X\\ . When collective effects are taken into account, the "bare" potential of interaction Ak.k'(J,J') 
is replaced by the "dressed" potential of interaction l/£'k,k'(J, J', w), without changing the overall structure of the 
kinetic equation. This is similar to the case of homogeneous plasmas where the bare potential iiCk) in the Landau 
equation is replaced by the "dressed" potential -u(k)/|e(k, k • v)| in the Lenard-Balescu equation. 

D. The relaxation time of the system as a whole 

The kinetic equation (|56l) is valid at the order 1/iV so it describes the "collisional" evolution of the system on 
a timescale ~ Nto. This kinetic equation conserves the total mass M = J f{J,t)dJ and the total energy E = 
J f{3,t)H{3)dJ. It also monotonically increases the Boltzmann entropy S = — /(//m)(J,i) ln(//m)(J,i) dJ {H- 
theorem). The Boltzmann distribution /e(J) = Ae"'^™^'^"'' is a steady state of this kinetic equation. The derivation of 
these results can be found in [63, ll^- The collisional evolution of the system (at the order 1/A^) is due to a condition 
of resonance between distant orbits. The condition of resonance, encapsulated in the (5-function, corresponds to 
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k' • J7(J') — k • r2(J) with (k', J') 7^ (k, J). If there are resonances during the whole evolution, the system will relax 
towards the Boltzmann distribution (provided that this equilibrium state exists). In that case, the relaxation time is 

tR - NtD. (61) 



This is the case, in general, for spatially inhomogeneous systems in any dimension [106[. This is also the case for 
spatially homogeneous systems in d > 1 dimensions. In that case, the condition of resonance corresponds to k-v' = k-v 
with v' 7^ V (the self-interaction at v = v' does not produce transport since the term in parenthesis in Eq. (|57p 
vanishes identically). However, for spatially homogeneous one dimensional systems, there is no resonance and the 



Lenard-Balescu collision term vanishes [68|, |80[. The kinetic equation ([57)) reduces to df/dt — 0, so the distribution 
function does not evolve at all on a timescale ^ Nti). In that case, the relaxation time is larger than Nt]j. We 
therefore expect that 

t_R > Nto (ID homogeneous). (62) 

Since the relaxation process is due to more complex correlations, we have to develop the kinetic theory at higher 
orders (taking three-body, four-body,... correlation functions into account) to obtain the relaxation time. If the 
collision term does not vanish at the next order of the expansion in powers of 1/A^, the kinetic theory would imply 
a relaxation time of the order N'^tj). However, the problem could be more complicated and yield a larger relaxation 
time (see the Conclusion for a more detailed discussion). In fact, it is not even granted that the system will ever 
relax towards statistical equilibrium. The evolution may be non-ergodic and the mixing by "collisions" inefhcient. 
By contrast, for spatially homogeneous systems in d > 1 dimensions [see Eq. (j57p ]. and for spatially inhomogeneous 
systems in any dimension [see Eq. (|56p]. since there are potentially more resonances, the relaxation time could be 
reduced and approach the natural scaling Nto predicted by the first order kinetic theory [63,|6g. However, very little 
is known concerning the properties of Eq. (|56p and its convergence (or not) towards the Boltzmann distribution. It 
could approach the Boltzmann distribution (since entropy increases), without reaching it exactly if the resonances 
disappear at some point of the evolution (see Refs. [63, [g^ for more details). 

III. STOCHASTIC PROCESS OF A TEST PARTICLE: DIFFUSION AND FRICTION 

A. The Fokker-Planck equation 

In the previous section, we have studied the evolution of the system as a whole. We now consider the relaxation 
of a test particle is a bath of field particles with a steady distribution /(J). We assume that the field particles are 
either at statistical equilibrium with the Boltzmann distribution (thermal bath) , in which case their distribution does 
not change at all, or in a stable steady state of the Vlasov equation with a profile that forbids any resonance (as we 
have just seen, this profile does not change on a timescale of order Ntn). We assume that the test particle has an 
initial action Jq and we study how it progressively acquires the distribution of the bath due to "collisions" with the 
field particles. As we shall see, the test particle has a stochastic motion and the evolution of the distribution function 
P(J, i), the probability density that the test particle has an action J at time t, is governed by a Fokker-Planck equation 
involving a diffusion term and a friction term that can be analytically obtained. The Fokker-Planck equation may 
then be solved with the initial condition P(J,0) = (5(J — Jq) to yield P{3,t). This problem has been investigated in 
our previous papers [67l . l68j , but we shall give here a direct and more rigorous derivation of the coefficients of diffusion 
and friction, taking collective effects into account. 

The equations of motion of the test particle, written with angle-action variables, are 

dw _, , 9(5$, , , dJ dS<P, , , 

-^n(J) + — (w,J,t), -^- — (w,J,0. (63) 

They include the effect of the mean field which produces a zeroth order motion characterized by the pulsation J7(J) 
plus a stochastic component S^ of order 1/\/N (fluctuations) which takes the deviations from the mean field into 
account. They can be formally integrated into 

m{t)=^ + j n{J{t'))dt' + J ^(w(i'),J(t'),Orfi', (64) 

/■* (9(5$ 
J(t)=J- / ^(w(i'),J(i'),0< (65) 

Jq (7W 
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where we have assumed that, initially, the test particle is at (w, J). Note that the "initial" time considered here does 
not necessarily coincide with the original time mentioned above. Since the fluctuation (5$ of the potential is a small 
quantity, the foregoing equations can be solved iteratively. At the order 1/N , which corresponds to quadratic order 
in (5$, we get for the action 



Mt) = J, - 



[ dt'^{yv + m\j,t') + ^ [ dt' I dt" [ di"'^^(w + ot',j,t')^(w + m"',j,r) 

Jo owi dJj Jo Jq Jo dwidwi dwj 

+ dt' dt" " (w + »t',J,tO^(w + »t^J,t") 
Jo Jo dwidJj dwj 

- dt' dt",f^(w + m',J,t')^(w + m",J,i")- (66) 

Jo Jo dwidw-j dJj 

As the changes in the action are small, the dynamics of the test particle can be represented by a stochastic process 
governed by a Fokker-Planck equation [8l|. If we denote by P(J, t) the probability density that the test particle has 
an action J at time t, the general form of this equation is 

dP 8^ 8 

The diffusion tensor and the friction force are defined by 



1 



A,(J) - ^ lim -HMt) - J.,)[J,{t) - J,)), (68) 



A,(J)= lim i(J,(i)_J,). (69) 

t-> + oo t 

In writing these limits, we have implicitly assumed that the time t is long with respect to the fluctuation time but 
short with respect to the relaxation time (of order Ntjj), so the expression (|55)) can be used to evaluate Eqs. (p5)) and 
(|69l) . As shown in our previous papers [63, IsM] , it is relevant to rewrite the Fokker-Planck equation in the alternative 
form 

The total friction is 

dn 

A^^pvoi^au^^ (71) 

where Fpo; is the friction due to the polarization, while the second term is due to the variations of the diffusion 
tensor with J. As we shall see, this decomposition arises naturally in the following analysis. The two expressions 
(p7)) and ([70)1 have their own interest. The expression (p7)) where the diffusion tensor is placed after the second 
derivative 8^ {DP) involves the total friction A and the expression ([70)1 where the diffusion tensor is placed between 
the derivatives dDdP isolates the friction by polarization Fpoi- We shall see in Sec. IIII El that this second form is 
directly related to the Lenard-Balescu-type equation ([56| . It has therefore a clear physical interpretation. We now 
calculate the diffusion tensor and the friction force from Eqs. (|68p and ([M|) . using the results of Sec. Ill CI that allow 
to take collective effects into account. The reader not interested in the details of the calculations may directly go to 
the final results summarized in Sec. IIII El 

B. The diffusion tensor 

We first compute the diffusion tensor defined by Eq. (p5|) . Using Eq. ([M|) . we see that it is given, at the order 
1/iV, by 



1 /■* /■* / 86^ 8S^ 

A, = :^ / ^^' / dt" (^(w + f2<',J,t')^(w + f2<",J,0 
2iJo Jo \8w, dwj 



(72) 
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By the inverse Fourier-Laplace transform, we have 

dS^ , „., , .,,oo-i>, „.„ , .„, \ ■^r^ I auj ^ir^ i aui , , , ,k-(w+f2t')g-i"t'gik'-(w+ot")g-icj't" 



a.,-'---' ^''''If <- + -"'^''">) -E/jEi^^.^" 



x(J$(k,J,a;)jl>(k',J,a;')). (73) 



Since the diffusion tensor depends only on the action, it can be averaged over w without loss of information. This 
brings a Kronecker factor Jk.-k' which amounts to taking k' = — k. Then, substituting Eq. P5|) in Eq. ([75]) . and 
carrying out the integrals over uj' and a;, we end up with the result 

(^'-^"''^■'''g'- + "'''-'''"')-'^''^"\Sy''^'^-'^^--'^"-""''''-'" V.,.(J.j'.i.--n')P ^'^''- 

(74) 

This expression shows that the correlation function appearing in Eq. (j72p is an even function of t' — t" . Using the 
identity 

dt' f dt" fit' - t") ^2 f dt' f dt" fit' - t") = 2 f dsit- s)fis), (75) 

"'0 "'0 "'0 Jo 



we find, for t — > +oo, that 



/■+°° /dS<P dS<P \ 

D.,^l ( — (w,J,0) — (w + n.,J,.s))ds. (76) 



(77) 



This is the Kubo formula for our problem. Replacing the correlation function by its expression ()74p . we get 

/■+00 f 1 

A. = i2.rml ds gy .J'..fc,e'(— - ''jU^M^irjTTTrWfin 



Making the change of variables s ^ —s and (k, k') -^ (— k, — k'), we see that we can replace /„ ds by (1/2) /_^ ds. 
Then, using the identity 



6ico)^ I e-*|^, (78) 



— OO 



we obtain the final expression 

A, - ^(2^)^m ^ J dJ' k,k,j^—^j-j^-^^-^Sik . O - k' . 0')./(J'). (79) 

C. The friction due to the polarization 

We now compute the friction term defined by Eq. (|69| . We need to keep terms up to order 1/N. From Eq. (|66p . 
the first term to compute is 



1 /■* /aj£, 

1 ,/o \ 9w 



Ai = -- dt'(^—iv, + nt',J,t')). (80) 



By the inverse Fourier-Laplace transform, we have 



dw 

Using Eq. dST]), we find that 



^'^*(w + nt',J,t')) =«^ /f^ke"'-^"'+''*'^e-'"*'(^|.(k,J,c^)). (81) 



(<5l.(k,J,c.))=-(27r)''V /dJ' ,1 ^ ^¥^^y''"^! ■ (82) 

^ ^ ' ^^ ^ ^ 4:^7 Ak'(J,J',c^)i(k'-0'-w) ^ ^ 
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Now, using the fact that the test particle is initially located in (w, J), so that {Sf{w', J', 0)) = m5{w' — w)(5(J' — J), 
we obtain from Eq. ([20]) the result 



{Sf{k',3',0)) = j^me-^^--6{J'-3). (83) 

Substituting these expressions in Eq. (|8ip . and averaging over w (which amounts to replacing k' by k), we get 



k 

Therefore, the fr-iction term (|80p is given by 



We now use the Landau prescription w ^ w + iO+ and the Plemelj formula (|46|) to evaluate the integral over oj. The 
term corresponding to the imaginary part in the Plemelj formula is 






(86) 



Integrating over u and i', we obtain 



<'-?i:""(iwim)' <«'> 



The term corresponding to the real part in the Plemelj formula is 






' tJo ^ J^oo 27r 



/o — J-oo ^" Dk,k(J, J,w) k- $7- w' 

Integrating over i', we can convert this expression to the form 

/-|-oo J 1 11 
^^ n at wv o TiT{*sin[(k-f7-a.)t]+cos[(k.f2-^)t]-l}. (89) 

_oo 27r Z;k,k(J, J,w) (k- Si- w)^ t 

For t — > +00, using the identity 



lim "Zl = nSix), (90) 



1 — cos(te) 

t->+oo tx 

and integrating over w, we find that A} is given by Eq. (|87)) . just like A^° . Therefore, writing A/ = A^" + Aj = 
2Aj , we obtain 

A/ = -mVkIm[ ,^\ , A . (91) 

^ V^k,k(J,J,k-0)y ^ ' 

Finally, using Eq. ([M)) . we find that 

A/ = 7r(27r)''m ^ /" dJ' k ( k' • ^ j 5(k • fi - k' • fi' 



kk.^ V -/ ,i?k,k'(J,J',k.f2)|2- 



(92) 



As we shall see. A/ corresponds to the friction due to the polarization denoted Fpo; in Eq. ([7T|) . 

Remark: We can obtain Eq. ([M)) in a slightly more direct manner from Eq. (|84p by using the contour of integration 
shown in Fig. 9 of [73 • In that case, the integral over uj is just — 27ri times the sum of the residues at the poles of the 
integrand in Eq. (|84p . The poles corresponding to the zeros of Z3k,k(J, J, w) give a contribution that rapidly decays 
with time since Im(a;) < (the system is Vlasov stable). Keeping only the contribution of the pole a; = k • Jl that 
does not decay in time we obtain 

^{w + nt',J,t')\ = -imY^ , \ , ^, ■ (93) 

dw"- ' ' V 4- Z?k,k(J,J,k-0) ^ ' 



k 



Substituting this result in Eq. (|5(7)) . we get Eq. ([?T|) then Eq. 
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D. The part of the friction due to the inhomogeneity of the diffusion coefficient 

In the evaluation of the total friction, at the order 1/A^, the second term to compute is 



t Jo Jo Jo 9Jj \dwidwi ' ' dwj 

By the inverse Fourier-Laplace transform, we have 



.(w+f,*'.a,o|f(w+f.r,j,n) = -.Ei£Ei^^M 



e ^ 'e 



dwidwi 

xe^'^'(^+"*"')e-'"'*"'((5l'(k,J,w)jl>(k',J,w')). (95) 

Substituting Eq. (|43| in Eq. (|95|) . averaging over w (which amounts to replacing k' by — k), and carrying out the 
integrals over w' and w, wc end up with the result 

d S^ , „., , .,,9(5$, „.,„ , .,„A . ,„ ,w \r^ /■ .,,, , , i(kf2-k'r2')(t'-t"') 



-(w + Oi',J,i')?^(w + m"',J,i"')\ =i(27r)''mV f dJ'k,k,kie' 



dwtdwi , ^^^ 



x^--k7-7^/(J'). (96) 



dt" / dt'" f{t' - t'") = / dt'" it' - t"')f{t' - t'"), (97) 

Jo Jo 



|i?k,k'(J,J',k'-n')P" 

This expression shows that the correlation function appearing in Eq. (|94p is an odd function of t' — t'" . Using the 
identity 

Jo 
we find that 

,4." ^ . i2.fj- [ ^ ll <«"' g /.r .,.,.,g«' - '"'K"-"-^-"^"-'^'' |^_^_^,(,;,,.„,,|. /(J-). (98) 
This can be rewritten 

or, equivalently. 



k,k' 



,2.)Vl/'../* .."'g/..'y,.*— ..-.-. ^ ^______ I /p',. (100, 



Since the integrand only depends on i' — t'" ^ using the identity (|75p . we obtain for t — >■ +oo, 

r' + OO 



-(2.,-,„f ~.,g/^'..V<— "•..A (__L__) „,-, 



(101) 



Making the change of variables s -^ ~s and (k, k') -^ (— k, — k'), we see that we can replace /„ ds by (1/2) j_^ ds. 
Then, using the identity (|78p . we obtain the expression 

A- = .(2.)'^™ A g j ,, fc.fc.- |^^^^,(j^;,,,.^,)|, ^(k . o - k' . o')/(j') 



In the first term, we recover the diffusion coefficient (1791) so that 



J J _ U^tj 



'^^^J _/o_Mi_V^ fJ^lu^.s:t^, n, iJ i~,/\ ^ 



Tr{2n fmY^ f dJ' hkjd{k ■ U - k' ■ Ct'} 



djj "• ' ^J ■> ^ 'ajj VI-Dk,k'(J,J',k'-r2' 
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/(J'). (103) 



Finally, the third and fourth terms to compute are 

rt ft' I p,2 



'^"^iL "'l '«"(^<"+"''-'-'''l?<"' +"'"•■'•' 



(104) 



and 



Substituting the inverse Fourier-Laplace transform of the fluctuations of the potential in these equations, averaging 
over w (which amounts to replacing k' by — k) and summing the resulting expressions, we obtain 

Af'+^r = ] fdt' f df'Y: I P / ^fc^fc^e'^^^^'e— *'e-'^"*"e-'"'*" A(5$(k,j,^)j|.(_k,j,^')). (i06) 



Substituting Eq. (|43p in Eq. (|lQ6p . and carrying out the integrals over cj' and w, we end up with the resuh 

This expression shows that the correlation function appearing under the integral sign only depends on t' — t" . Using 
the identity ((75|) we find, for t -^ +oo. that 

A;" + A'- = (2.)-..f".J.g/^J'^.>^.<."-"--'-^( |^^^^(,,;,,.„,)|. )/W. (108) 

Making the change of variables s ^ —s and (k, k') -^ (— k, — k'), we see that we can replace J^ ds by (1/2) J_^ ds. 
Then, using the identity (|78p. we obtain the expression 



Al- + A- = <2.rm g / dJ' k.k, A ( |^^^^(j;,,,.^,)|, ) .(k . n - k' . n')/(J'). (109) 

Finally, summing Eqs. (|103p and (|109p . we get 

Aii + Ajjj + Ajv = ^. (110) 

dJj 

E. Connection between the kinetic equation (|56|) and the Fokker-Planck equation (|70|) 

We have established that the diffusion tensor and the friction force are given by 

A, = ^(2^)^m ^ J dJ' hk, ^^^ ^^^j^ j,^ ^, ^,^^^ S{k . O - k' . 0')/(J'). 



(Ill) 



k,k' 



1 ... .. ,j ..,.(., 9f\ , dD 



A. = .(2.)^^^/.J-fc. ^^^^^^^jj,^^^^^, .(k.»-k-.O0(k- 



'u 



aj' / a J, 



(112) 



Comparing Eq. (|112p with Eq. ([7T]) . we see that the friction by polarization is 
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On the other hand, using an integration by part in the first term of Eq. (|112p . the total friction can be written 

Finally, using Eqs. (|llip and (|113p . we find that the Fokkcr-Planck equation ([7(I| becomes 



k.k' 



When collective effects are neglected, we recover the results obtained in [671 l68l | by a different method (related 
expressions of the diffusion and friction terms have also been obtained in [59l466| ). For spatially homogeneous systems. 



we recover the results of Hubbard [ll| (see also [37| when collective effects are neglected). As observed in our previous 
works, we note that the form of Eq. (jllSp is very similar to the form of Eq. ([56| . This shows that the Fokker- 
Planck equation (jlisp . with the diffusion coefficient ()llip and the friction (|113p . can be directly obtained from the 
Lenard-Balescu-type equation (|56p by replacing the time dependent distribution /(J',i) by the static distribution 
/(J') of the bath. This procedure transforms an integro-differential equation (|56p into a differential equation (jlisp . 
Although natural, the rigorous justification of this procedure requires the detailed calculation that we have given here. 
In fact, we can understand this result in the following manner. Equations (|56p and (jllSp govern the evolution of 
the distribution function of a test particle (described by the coordinate J) interacting with field particles (described 
by the running coordinate J'). In Eq. ([56|. all the particles are equivalent so the distribution of the field particles 
f{3',t) changes with time exactly like the distribution of the test particle f{3,t). In Eq. (|115p . the test particle and 
the field particles are not equivalent since the field particles form a "bath" . The field particles have a steady (given) 
distribution /(J') while the distribution of the test particle f{J,t) = NmP{J,t) changes with time. 

F. Thermal bath: Boltzmann distribution 

For a thermal bath, the field particles have the Boltzmann distribution of statistical equilibrium 

/(J) = Ae"'3™^(-'\ (116) 

where H{J) is the individual energy of the particles. We note the identity 

^ = -pmf{3)n{3), (117) 

where we have used dH/d3 = Q,{3). Substituting this relation in Eq. (|113p . we obtain 

Fp„, = -7r(27r)'*™2^^|dJ'k^^— ^j^l-^— ^<5(k. n -k' • n')(k' • n')f{3'). (118) 

Using the (5- function to replace k' • Jl' by k • O, and comparing the resulting expression with Eq. (|llip . we find that 

Ff°' = -PmD,j{3)Vij{3). (119) 

This can be viewed as a generalized Einstein relation connecting the friction force to the diffusion tensor (fiuctuation- 
dissipation theorem). We stress that the Einstein relation is valid for the friction by polarization Fpoj, not for the 
total friction A that has a more complex expression due to the term djDij. We do not have this subtlety for the usual 
Brownian motion where the diffusion coefficient is constant. For a thermal bath, using Eq. (|119p . the Fokker-Planck 
equation ([70)1 takes the form 



dP d 

'dt ~'dli 



dP 

3 



D,j{3)[—+PmPnj{3) 



(120) 



where Dij{3) is given by Eq. (jllip with Eg. (|116|) . Recalling that ri(J) = dH/d3, this equation is similar to 
the Kramers equation in Brownian theory [l9J. This is a drift-diffusion equation describing the evolution of the 
distribution P{3,t) of the test particle in an "effective potential" Ueff{3) = H{3) produced by the field particles. For 
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t — ?► +CXD, the distribution of the test particle relaxes towards the Boltzmann distribution Pe(J) = {A/Nm)e~^"^^'^'^\ 
Since the Fokker-Planck is valid at the order 1/A^, the relaxation time scales like 

4"*'^ - NtD- (121) 

If the field particles are in a stable steady state of the Vlasov equation that forbids any resonance (see Sect. [ IIPp , 
the diffusion tensor and the friction force can be simplified in a manner similar to that described in |67l . |68|. In 
that case, the test particle still relaxes towards the distribution /(J) of the bath on a timescale i^*'' ^ Nto- This 
timescale is shorter than the relaxation time of the system as a whole tn > Nto so that the bath approximation, 
which assumes that the distribution of the field particles is "frozen" , is justified. 

IV. CONCLUSION 

In this paper, we have developed a kinetic theory of systems with long-range interactions which takes collective 
effects and spatial inhomogeneity into account. This improves our previous works [63, IS^ where collective effects 
were neglected (leading to the Landau- type kinetic equation (j59p ) or taken into account in an unsatisfactory manner 
(making an ad hoc factorization assumption) . The use of a biorthogonal basis of pairs of density-potential functions 
solves the shortcomings of our previous derivation. We recovered and confirmed the Lenard-Balescu-type kinetic 
equation ([55t derived by Heyvaerts ^0| from the Liouville equation, using the BBGKY hierarchy truncated at the 
level of the two-body distribution function (i.e. neglecting three-body correlations). Our der ivati on, starting from 
the Klimontovich equation and using a quasilinear approximation, may be technically simpler [l07l |. Furthermore, it 
is completely parallel to the one exposed in [2, U^ UM to derive the ordinary Lenard-Balescu equation for spatially 
homogeneous systems and in |5ll . |53| to derive a Lenard-Balescu-type equation appropriate to two-dimensional point 
vortices with axisymmetric distribution. Our approach therefore offers a unified description of kinetic theories for 
systems with long-range interactions. Interestingly, the same formalism can be used to describe the relaxation of a 
test particle in a bath of field particles through a Fokker-Planck equation. This is another advantage of our approach. 
We have obtained the expressions of the diffusion tensor and friction force directly from the equations of motion with 
angle-action variables, taking collective effects into account. We have shown that they can also be obtained from the 
Lenard-Balescu-type kinetic equation (j56p by making a bath approximation leading to the Fokker-Planck equation 
(|115p . For simplicity of notation, we have considered a single species system, but the general case of a multi-species 
system can be treated easily and leads to the same results as those obtained by Heyvaerts |69| . 

The kinetic theory developed in this paper is valid at the order 1/A^ so it describes the evolution of the system 
on a timescale Ntj:,. The dynamical evolution of the system is due to a condition of resonance between distant 
orbits. For spatially homogeneous one dimensional systems, there is no resonance so the Lenard-Balescu collision 
term vanishes (this result is well-known in plasma physics [49|, |50[ and it has been rediscovered recently for the HMF 
model [43, |4^). This implies that the relaxation time is larger than Nto- In that case, we need to develop the kinetic 
theory at higher orders (taking three-body, four-body,... correlations into account). At present, no such theory 
exists and the scaling of the relaxation time with N remains an open problem. We can try to make speculations 
based on available results of numerical simulations. The most natural scaling would be N'^tjj which corresponds to 
the next order term in the expansion of the basic equations of the kinetic theory in powers of 1/N |68| . An N'^ 
scaling is indeed observed numerically |8J, |85| for spatially homogeneous one-dimensional plasmas (in that context, 
N represents the number of charges in the Debye sphere). However, the scaling of the relaxation time may be more 
complex. For example, for the permanently spatially homogeneous HMF model, Campa et al. [8g report a relaxation 
time scaling like e^tn- This timescale is, however, questioned in recent works [83| who find a N'^to scaling. This 
would be more natural on a theoretical point of view. It could also happen that the system never reaches statistical 
equilibrium, i.e. the evolution may be non-ergodic and the mixing by "collisions" inefficient. This may be the case 
for the two-dimensional point vortex gas when the profile of angular velocity is monotonic |5ll453l |. In that case, 
there is no resonance and the relaxation time is larger than Nto- However, the precise scaling with N is still an 
open problem, and even the relaxation towards the Boltzmann distribution is uncertain. The situation is different for 
spatially homogeneous systems in d > 1 dimensions and for spatially inhomogeneous systems in any dimension (and 
for non-axisymmetric configurations of the point vortex gas). In these cases, there are potentially more resonances 
so the relaxation time can be reduced and achieve the natural scaling Nto corresponding to the first order of the 
kinetic theorv | 68 | . An Nto scaling is indeed observed numerically for spatially inhomogeneous one dimensional stellar 
systems [6^, |88| - [9(]| and for the spatially inhomogeneous HMF model [91| . On the other hand, for the HMF model, if 
an initially spatially homogeneous distribution function becomes Vlasov unstable during the collisional evolution, a 
dynamical phase transition from a non- magnetized to a magnetized state takes place (as theoretically studied in [92|) 
and the relaxation time could be intermediate between N'^to (permanently homogeneous) and Nto (permanently 
inhomogeneous). In that situation, Yamaguchi et al. [93] find a relaxation time scaling like N^tjj with 5 = 1.7. 
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The previous argument (suggesting 1 < S < 2) may provide a first step towards the explanation of tliis anomalous 
exponent. 

The results of this paper also apply to self-gravitating systems that are, by essence, spatially inhomogeneous and 
limited in extension. The consideration of spatial inhomogeneity avoids the logarithmic divergence occurring at large 
scales when one makes the local approximation [70|. There remains, however, a logarithmic divergence at small scales 
due to the neglect of strong collisions in the kinetic theory. On the other hand, self-gravitating systems are very 
particular since no statistical equilibrium state exists in a strict sense (the entropy is not bounded from above) |3|-[5| . 
Therefore, the Lenard-Balescu-type kinetic equation (|56p does not relax towards a steady state. On the contrary, the 
system takes a "core-halo" structure and keeps evolving. The halo expands due to the evaporation of high energy stars 
j3G | and, when the system becomes sufficiently centrally condensed, the core collapses as a result of the gravothermal 
catastrophe [3l|, [32| caused by the negative specific heat of the central region. This evolution is clearly illustrated in 



the numerical simulation of Cohn [28 
leads to the formation of binary stars 



based on the orbit-averaged Fokker-Planck equation. Core collapse ultim ately 
27| , followed by a post-collapse evolution [9J] , and gravothermal oscillations [93 . 
The formation of binary stars is not taken into account in the kinetic theory because it results from triple collisions 
(three-body correlations) and is of the strong interaction type. We also note that, since the Boltzmann distribution is 
not normalizable for self-gravitating systems (infinite mass problem), the thermal bath approximation that we have 
developed in Sec. IIIIFl is not justified. However, it could be easily generalized by replacing the Boltzmann distribution 
by the King model which takes evaporation into account and leads to a truncated Boltzmann distribution with a finite 
mass viewed as a metastable equilibrium state [5| . 

At the physical level, one interest of the kinetic theory is to show that the "collisional" evolution of the system is 
controlled by resonances encapsulated in a (5-function. This explains the scaling of the relaxation time with N. For 
example, in the case of one dimensional potentials, the relaxation time of spatially homogeneous systems is larger than 
Ntr) because there is no resonance, while it is of order Ntu for spatially inhomogeneous systems due to the arising of 
resonances. Obtaining this result was our original motivation to develop a kinetic theory for spatially inhomogeneous 
systems in [63, [G^ . At a practical level, we believe that the Lenard-Balescu equation (|55|) with angle-action variables 
is very complicated to solve, so its practical interest is limited. With present-day computers, it is much easier to solve 
the A^-body problem directly. Although the description of the motion becomes simple in angle-action variables, all the 
complexity is now embodied in the relation between position-momentum variables and angle-action variables. The 
usefulness of equation (j56|) is therefore limited to systems for which this relation can be established, either analytically 
or numerically. Some simplification can be gained if we restrict ourselves to spherically symmetric stellar systems as 
done by Heyvaerts [631 for stellar systems. But, even in that case the final kinetic equations remain complicated and 
approximations may be welcome to simplify the problem and lead to more tractable equations. The Fokker-Planck 
treatment that we have given here is one possible approximation. The restriction to one dimensional systems with 
simple potentials (e.g. the HMF model) is another one. 

Appendix A: Bare and dressed potentials of interaction 

We assume that the particles interact via a long-range binary potential w(|r — r'|) so that 

$(r, i) = y u(|r - r'|)/(r', v', t) dr'dV. (Al) 

We introduce angle-action variables and decompose the potential of interaction in Fourier modes 

u{\r{w, J) - r'(w', J')|) = Y. e'^^^-^-'^'-^') Ak,k' (J, J'), (A2) 

k.k' 

where the Fourier coefficients are given by 

AuM' (J, J') = -j^ J «(|r(w, J) - r'(w', JODe-'^-^— '^'-') dwdw' . (A3) 

The fluctuations of the potential are related to the fluctuations of the distribution function by 

S<P{w, J, i) = f u{\r - r'|)(5/(w', J', t) dw'dJ', (A4) 

where we have used the identity dr'dv' = dw'dJ' valid for a canonical change of variables. Taking the Fourier-Laplace 
transform of this expression, and using Eqs. (fTSj) and (jA3p . we get 



6^{k, J, cj) = (277)"^ J2 I ^J'^k,k' (J, 3')5f{\i', J', ijj). (A5) 
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Substituting Eq. (1211) in Eq. (|A5|) . wc obtain the Fredholm integral equation 

J df 






k' " k 



This equation relates the Fourier-Laplace transform of the fluctuations of the potential to the Fourier transform of 
the initial fluctuations of the distribution function. It is therefore equivalent to Eq. (|37| . When collective effects are 
neglected, the foregoing expression reduces to 

<5$(k, J,a.) = (2^)'^^|dJ'Ak,k'(J, J')^^^^^^^- (A7) 

Comparing this expression with Eq. ([57| . we see that —1/D]f^_]f^r{3,3' ,uj)bare is just the Fourier transform Ak,k'(J, J') 
of the "bare" potential of interaction u with angle-action variables. When collective effects are taken into account, 
Eq. (|A7P is replaced by Eq. ([57)1 where according to (|A6|) the "dressed" potential of interaction satisfies 



.11 df 



(2ttY'Y [ d3"Akk„(J,J") ^"'93" i ^ = -Akk'(J,J')- 

f^i ^k"-n"-wDk"k' J",J',w) ■ 



(A8) 



£'k,k'(J,J','^) i^J ^k"-n"-c^£)k".k'(J",J',t^) 

Appendix B: The resolvant 

Substituting Eq. ^7^ in Eq. ([21]), we obtain 

This equation relates the Fourier-Laplace transform of the fluctuations of the distribution function to the Fourier 
transform of the initial fluctuations of the distribution function. It can be written 

<5/(k, J,..) = Y. /rfJ"i?k,k"(J, J",^)^", J",0), (B2) 

k" -^ 

where 

i?k.k"(J, J",c.) = - Ji-i_(2^)^^ ^ 1 ^ + JJ£_i!) 4_,„, (B3) 

is called the resolvant. If we consider an initial condition of the form 5/(w, J, 0) ~ m6(\v — w')(5(J — J') implying 

5/(k, J,0) = ^e-k.w'^(j _ j,^^ (B4^ 

we find that 

k • — rr I ^ 1 777 / 1 

This quantity is called the propagator of the linearized Vlasov equation. If we substitute Eq. (|A5[) in Eq. (PT|) , we 
obtain the integral equation 

df{k,J,c.) = Ji-IL(2,r)'^^ /dJ'Ak,k'(J, J')W', J',^) + .ff^y'°\ . (B6) 

k-s2 — w ■^ J i(k-S2 — cj) 

k' 

This equation relates the Fourier-Laplace transform of the fluctuations of the distribution function to the Fourier 
transform of the initial fluctuations of the distribution function, so it is equivalent to Eq. (jBll) . 
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Appendix C: Auto-correlation of the fluctuations of the one- particle distribution 

According to Eq. ([20]) . we have 

(<5/(k, J, 0)6 fik', J', 0)) = J -^-^^e-^i^-+^---) {Sfiv,, J, 0)^/(w', J', 0)) 

^^^^^-.(k.w+k'.w') j^^^(^^ j^ o)^^(^,^ j^ 0)^ _ /(j);(j')] . (ci) 



The expression of the discrete distribution function recalled in the first paragraph of Sec. Ill Al leads to 

(/d(w,J,0)/d(w',J',0)) = m2^(5(w-w,)<5(J-J,)<5(w'-w,)<5(J'-Jj)) 



1-3 
2 



Y^ {S{w - Wi)S{J - Ji)'5(w - w')S(3 - J')) 



m 



Y^ ((5(w - Wi)S{J - J,)(5(w' - Wj)5(J' - Jj)) 

= m/(J)5(w-w')<5(J-J') + /(J)/(J'), (C2) 

where we have assumed that there is no correlation initially (if there are initial correlations, it can be shown that they 
are washed out rapidly [72| so they have no effect on the final form of the collision term). Combining Eqs. (JC1|) and 
([C2|) . we obtain 

(J/(k, J, 0)6 f{k', J', 0)) = ^4,^k'<5(J - J')mf{J). (C3) 
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